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We experimentally study a two-dimensional (2D) granular pack of photoelastic disks subject to vertical manual tap¬ 
ping. Using bright and dark field images, we employ gradient-based image analysis methods to analyse various structural 
quantities. These include the packing fraction (</>), force per disk (Fj) and force chain segment length (/) as a function 
of the tapping number (r). The increase in packing fraction with tapping number is found to exponentially approach 
an asymptotic value. An exponential distribution is observed for both the cumulative numbers of force per disk Fj : 
Ncum(Fd) - exp(-|4), and force chain segment length / : N cllm (l) = A ; exp(-2-). Whereas the coefficient A F varies 

with r for force per disk, the force chain segment length shows no dependence on r. The r dependence of Fj and (j> 
allows us to posit a linear relation between the total force of the granular pack (F* of (r)) and </>(r). 


1. Introduction 

Granular packing is ubiquitous in everyday life. It is com¬ 
mon knowledge that a denser granular pack can be achieved 
by tapping the pack. Clogged granular flow can be unjammed 
and structural foundations of buildings strengthened with tap¬ 
ping. Indeed, the first thing one does when in trouble with 
handling granular materials is to tap the container. Neverthe¬ 
less, the physical mechanisms concerning the effect of tapping 
on granular packs are not yet completely understood. Recent 
investigations on granular compaction have used the dimen¬ 
sionless maximum acceleration T = a max /g to characterize 
the strength of tapping and/or vibration applied to a granular 
pack, where a max is the maximum acceleration and g = 9.8 
m/s 2 is the gravitational acceleration. 1-11 Most previous stud¬ 
ies have used steady vibration to cause granular compaction. 
The final state attained by steady vibration is solely deter¬ 
mined by F value. 2 ' 3 ' 61011 The most efficient compaction is 
induced at a value of T ^ 2. 2 ' 3,7 When T is too small, the 
compaction takes long and grows logarithmically with time. 
When T is too large, on the other hand, it is difficult to at¬ 
tain the highly compacted state as a large amount of kinetic 
energy is delivered to the granular pack in such a strong vibra¬ 
tion. 7 However, granular compaction also depends on vibra¬ 
tion history. 4,5 Although steady vibration has a well-defined 
maximum acceleration value, it represents one particular in¬ 
stance of granular compaction. In general, natural vibration 
or tapping applied to granular pack is somewhat irregular. 
Hence, granular compaction induced by irregular perturba¬ 
tions such as manual tapping must be examined to understand 
compaction processes that occur in nature. 

To diagnose the physical mechanism of granular com¬ 
paction, an access to the inner stress structure created by gran¬ 
ular pack is necessary. In general, the granular packs exhibit 
inhomogeneous stress distribution which can be character¬ 


ized by a network of force chains. This force chain struc¬ 
ture is peculiar to granular assemblies and causes their com¬ 
plex rheological behaviors. The force chain structure can be 
visualized in two-dimensional (2D) case. Using a 2D pack 
of photoelastic disks, the force chains can be observed via 
the retardation due to stress-induced birefringence of photoe¬ 
lastic material. 12-18 Using photoelastic disks, force applied 
to each disk can be measured. 12,13 More recently, the ap¬ 
plied force has been decomposed into the normal and tan¬ 
gential components using a computational image matching 
method. 1416 Relations among the shearing, isotropic com¬ 
pression, and jamming have been intensively studied using 
photoelastic disks. 14-18 To the best of our knowledge, how¬ 
ever, the tapping-induced granular compaction has not been 
studied with photoelastic disks. Therefore, we carry out the 
experiment with photoelastic disks to approach the physics of 
granular compaction via manual tapping. 

In this paper, we report the details of experimental method¬ 
ology developed to study the granular compaction. A 2D 
granular pack consisting of photoelastic disks is constructed. 
Then, the manual taps are added to the granular pack. The 
evolution of packing fraction and force chain structure in the 
granular pack are characterized by image analysis of photoe¬ 
lastic disks. Using the analysed results, the relation between 
the packing fraction and force chain structure is discussed to 
reveal what happens in the compacted granular pack. 

2. Experiment 

2 .1 Setup 

The experimental setup as shown in Fig. 1, consists of 
a 2D experimental chamber constructed with acrylic plates 
along the front and back and held together by aluminium 
bars along the sides and bottom. The chamber dimensions 
are 0.3 x 0.5 x 0.011 m in height, width and thickness, 
respectively. An accelerometer (EMIC Corp. Model: 710-C) 
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Fig. 1. Top view of the optical setup of the experiment. The experiment is 
carried out in a dark room to prevent stray light. The distance between the 
light source and the camera is about 1 m to make uniform angles of incident 
light into camera. The 2D experimental chamber is put vertically in front of 
the light source which is attached to a circular polarizer. A snapshot of the 
chamber is taken with a CCD camera (Nikon D70). The circular polarizer in 
front of the camera is set at 90° (cross-polarisation mode) relative to the other. 
Two types of images are obtained with and without the circular polarizer in 
front of the camera. 


is mounted on the top right corner of the chamber to measure 
the maximum acceleration ( a max ) experienced during the 
experiment. The chamber is filled with a bidisperse (to avoid 
crystallisation) set of 350 large (diameter is 0.015 m) and 
700 small (diameter is 0.01 m) photoelastic disks of 0.01 
m thickness (Vishay Micro Measurements, PSM-4). The 
chamber is vertically placed between a circular polarised 
LED light source and a CCD camera (Nikon D70), which 
acquires two types of images 2000 x 3008 pixels in size, 
corresponding to a spatial resolution of 1.76 x 10 4 m/pixel 
(MPP). The camera is placed 1 m in front of the experimental 
chamber. First, a bright field image (Fig. 2a) of the granular 
pack is acquired to measure the packing fraction and the 
disk centres and diameters for estimation of force per disk. 
A second dark field image (Fig. 2b) is acquired by placing a 
second circular polariser between the experimental chamber 
and the camera in cross-polarisation mode. This image 
provides the photoelastic intensities of the granular force 
chains. Images are acquired under dark room conditions to 
minimise ambient noise from extraneous illumination. 


2.2 Experimental Protocol 

Prior to start of the experiment, an initial configuration of 
low packing fraction is generated. It is preferable that the ini¬ 
tial packing fraction be small since this study focuses on gran¬ 
ular compaction via manual tapping. However, when disks 
are introduced in a vertically standing chamber, initial com¬ 
paction occurs from disk impacts. Therefore, the disks were 



Fig. 2. (a) An example of bright field image by which the packing fraction 

and the position of photoelastic disks are obtained, (b) The corresponding 
dark field image by which the structure of the force chains is analysed 


introduced by spreading them in the chamber while it is hori¬ 
zontally laid down, and then the chamber was vertically fixed, 
thus assuring a small initial packing fraction. A pair of bright 
and dark field images is then acquired for this initial configu¬ 
ration. 

The system is then perturbed by providing a manual tap¬ 
ping to the experimental chamber. In particular, each manual 
tapping is defined as adding two impulses to each bottom edge 
of the experimental chamber. Whereas this tapping protocol is 
not systematically controlled as in the case of an electromag¬ 
netic shaker, for instance, it was specifically chosen to mimic 
the situation of stochastic impulse forcing observed in many 
natural processes. In any event, the accelerometer attached to 
the experimental chamber measures the acceleration experi¬ 
enced during tappings, from which dimensionless accelera¬ 
tion is defined as T = yy. The experiments reported here are 
in the regime of F -3-4. This tapping acceleration is large 
enough to achieve the efficient compaction. Following each 
manual tap, a pair of bright and dark field images is acquired 
for subsequent analysis to determine the evolution in the pack¬ 
ing fraction ((f)), force per disk (P,i), and the force chain seg¬ 
ment length (/) as a function of the tapping number r. Each 
experimental run consists of the initial configuration followed 
by nine manual taps (r = 9), thus providing ten pairs of bright 
and dark field images per run. Nine experimental runs under 
identical experimental conditions were conducted. 

2.3 Image Analysis 

Here we explain the image analysis methods employed to 
extract the packing fraction </>, force per disk E,j, and force 
chain segment length /, from the bright and dark field images. 
The image analysis software. Image J, 19 was used to analyse 
the experimental image data. 

2.3.1 Determination of Packing Fraction 

In this study, we define the packing fraction as f — 
S,/(S m +S V ), where S t is the theoretical total area of the pho¬ 
toelastic disks, S m is the total area of photoelastic disks mea- 
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Fig. 3. Disk identification and measurement of force per disk from image 
analysis. The area and centre of each disk are obtained from the bright field 
image (a). Sample disk center and circumference are shown in red for a large 
disk. This information is then used in the corresponding dark field image (b) 
to obtain the photoelastic signal at disk contacts, and subsequent analysis is 
employed to measure the force per disk. 


Fig. 4. Calibration method using a vertical ID chain of photoelastic disks. 
Bright (a) and dark (b) field images of the vertical ID chain. Image analysis 
of the ID chain images using algorithms explained in section 2.3.2 provided 
values of (G 2 ). Force per disk estimated from gravitational forcing F c was 
then used to relate F c versus (G 2 ) for both large (solid red circles) and small 
(solid blue triangles) disks (c). Solid lines through the calibration data are 
quadratic fits, which were used in experimental measurement of force per 
disk Fj. 


sured from the bright field images, and S ,, is the total void area 
measured from the bright field images. Whereas theoretically, 
S t = S m , in reality S m /S, - 1.1 due to the optical distortion 
of the image between the centre and edges of the bright field 
image (see Fig. 2a) and the thickness gap between disks and 
the chamber wall. When granular compaction occurs under 
manual tapping, whereas S m remains almost constant, S v de¬ 
creases due to reduction in area of the voids between disks. 
Therefore, a measurable increase in packing fraction is ob¬ 
served with tapping number r. Often, the packing fraction is 
calculated as the ratio of the area occupied by the photoelas¬ 
tic disks to the total chamber area. This definition is reason¬ 
able when the granular pack is enclosed on all sides. However, 
since this experiment is conducted with the upper side of the 
experimental chamber left open, an accurate estimation of the 
total chamber area is not possible. The same situation also 
arises in estimation of packing fraction for granular heaps or 
sand piles. 20 

For calculation of packing fraction the theoretical area 
of the disks S t was first calculated from the known number 
of large and small disks, whose diameters were already 
available, yielding S t = 1.17 x 10 1 m 2 for experiments 
reported here. For calculation of S m and S ,,, the bright field 
image was first binarized by using ImageJ which resulted 
in dark disks on a bright background. The pixel area of the 
dark regions multiplied by the spatial resolution (MPP) then 
provided the value of S m , and inversion of images was used 
to obtain the value of S v . 

2.3.2 Characterization of photoelastic intensity gradient 

Extant studies have used the photoelastic signal to mea¬ 
sure contact forces in one of two ways. The first method esti¬ 
mates force per disk by using the photoelastic intensity gradi¬ 
ent. 12-17 The second method estimates force per disk via com¬ 
putational image matching. 14 - 16 This study applies the former 


method for measurement of force per disk as the image reso¬ 
lution obtained is insufficient to measure forces by the latter 
computational matching scheme. The algorithm applied here 
for force measurement is similar to the one adapted by Howell 
et al.. 12 

For the given intensity I for each image pixel (8 bit, gray 
scale), Sobel filter was applied to obtain the squared hori¬ 
zontal (V//,) 2 and vertical (V/,,) 2 gradients in intensity. Their 
sum |V/| 2 = (V//,) 2 + (V/,,) 2 provides the squared gradient in 
intensity per pixel. The mean squared intensity gradient over 
all pixels within a disk was then defined as (G 2 ) = (|V/| 2 ). 
Computation of (G 2 ) on each disk requires knowledge of 
each disk centre and its area; information available from the 
bright field image (Fig. 3) which proceeded in three steps: 
(1) binarizing a bright field image, (2) splitting disk areas 
of contiguous binarized intensity into individual disks, and 
(3) measuring each disk centre position and area. Step 1 is 
identical to the packing fraction measurement method. In 
step 2, a watershed algorithm was employed to discriminate 
between sharp gradients in intensity among disks, usually 
referred to as mountain (low intensity gradient) and river 
(high intensity gradient), and separates out individual disks. 
This is necessary to identify the disk perimeters along which 
the photoelastic intensities of contact forces exist. Following 
this watershed procedure, each disk centre position and 
area were measured in step 3. Applying these results on 
the dark field image, the mean squared intensity gradient in 
photoelastic signal (G 2 ) was then obtained for each disk. 

2.3.3 Force Calibration 

The force per disk was calibrated using a vertical one¬ 
dimensional (ID) chain and the measurement method 
explained in section 2.3.2 to obtain force calibration curves 
which convert (G 2 ) to force. A vertical ID chain of photoe- 
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Original image Binarized Thinned 


Fig. 5. Method of stress chain thinning and definition of stress chain seg¬ 
ment length. A thinned stress chain (c) is obtained by binarizing original 
image ((a) —» ( b )), and thinning it ((b) — » (c)). A stress chain segment length 
is defined by a liner distance between intersections or end points in thinned 
stress chain image. 


lastic disks of height 0.3 m as shown in Fig. 4a & b was 
constructed. The ID chain consisted of either 20 large disks 
or 30 small disks. A pair of bright and dark field images was 
then obtained, and the image analysis methods (section 2.3.2) 
were applied to calculate (G 2 ) for each photoelastic disk. 
F c {n) in Newtons (where n is the position of the disk from 
the top of the ID chain), applied force per disk in the vertical 
ID chain was estimated from the relation F c = n x mg , 
where m is the mass per disk (1.8 x 10~ 3 kg for large disk 
and 0.8 x 10 3 kg for small disk, respectively). In Fig. 4b, 
fringes on the boundary between disks and sidewalls cannot 
be observed in the dark field image. Therefore, the effect 
of sidewalls was neglected in calibration measurements. 
Figure 4c shows the calibration data obtained for both disk 
sizes. The quadratic fits of the calibration data were then 
used as the final calibration curves for measurements of 
force per disk in the experimental data. Since the adopted 
procedure does not involve computational image matching of 
photoelastic fringes, only the total force applied to disk can 
be measured in this study, and cannot be decomposed into the 
normal and tangential components. 


2.3.4 Force chain segment length measurement 

The segment length l of force chains forms one of the 
structural variables measured in this experimental study. We 
employed a standard image analysis technique known as 
thinning method, an example of which is shown in Fig. 5. 
A dark field image (Fig. 5a) was binarized (Fig. 5b) and a 
skeletonize procedure (also known as the erosion method 
or the bleeding algorithm) in ImageJ was used to thin the 
segment down to a line of single pixel thickness. The force 
chain segment length was then defined as the linear distance 
between intersections or end points of the chain in the thinned 
force chain image (Fig. 5c). 



Fig. 6. The variation of packing fraction with manual tapping. Packing 
fraction increases with each manual tapping and approaches the steady state 
(0o = 0.79). Mean value of nine runs is reported, and the error bars rep¬ 
resent the standard error of nine runs. The dotted curve is the fit, (p(r) = 
00 + Aexp(-^), where 0o = 0.79, A - -1.39 x 10 -2 and to = 2.27 are 
values obtained for the fit parameters. 


3. Results 

3.1 Packing Fraction 

The calculated results for packing fraction at each tapping 
number r are shown in Fig. 6, where r = 0 represents the ini¬ 
tial configuration. The experimental data for (Mr) are fit with 
the function 0(r) = tpo +A exp(-^), where <f> 0 , A and r 0 are fit 
parameters. The fit parameter values for this study were found 
to be 0o = 0.79, A = -1.39xl0~ 2 and to = 2.27, respectively. 
As a result. Fig. 6 reveals that the packing fraction exponen¬ 
tially asymptotes towards a final steady state packing fraction 
value, in agreement with prior results reported by Bandi et 
al.. 15 Figure 6 reports the mean over nine experimental runs, 
with the error bars being their standard error. 

3.2 Force per disk 

The force on each disk in the granular pack was measured 
using the method described in sections 2.3.2 & 2.3.3. Figure 7 
shows the cumulative number distribution of force per disk 
at each tapping number r in the granular pack. The range of 
force per disk F d in Fig. 7 is wider than the calibration range 
(Fig. 4c). However, the calibration is performed under ID di¬ 
ametral compression, i.e., coordination number is 2 in the cal¬ 
ibration. In the granular pack, on the other hand, average co¬ 
ordination number is almost 4. Thus, the force applied to each 
disk can approximately be two times greater than the calibra¬ 
tion case. Then, the force magnitude of each contact point in 
the granular pack is almost within the calibration range.The 
distribution can be approximated by the exponential form: 

Ncum(F d ) = A/rexp^-^j, (1) 
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Fig. 7. Cumulative number distributions of force per disk at each tapping 
number r in log-linear scale. The black solid line represents the initial con¬ 
figuration (r = 0) whereas colored lines represent the respective compacted 
states for various r values. The data represent the mean value of nine experi¬ 
mental runs. 


Fig. 8. The variation of fit parameter Af as a function of tapping number 
r. Af increases with each manual tapping. The Af values are calculated from 
the fitting by Eq. (1) with fixed Fq ( = 0.09 [N]). The error bars represent the 
uncertainty of the fitting. 


where N cum (Fd) is the number of disks on which the applied 
force is equal to or greater than Fd. A f and Fq are fit 
parameters with Fq having dimensions of force. Figure 7 is 
obtained from the mean value of nine identical experimental 
runs, and exhibits a roughly exponential distribution with 
almost constant slopes across all values of r, initial as well 
as the final compact state. The fit parameters were found to 
be A f = 1.53 x 10 3 and Fq = 9.06 x 1CT 2 N at the initial 
configuration (r = 0). This result suggests that the functional 
form of the cumulative force distribution itself is invariant to 
the compaction under manual tapping as it yields the same 
slope for the exponential tail for all r values. This result 
is qualitatively consistent with the previous study in which 
Liu et al. and Coppersmith et al. measured the cumulative 
distribution of force exerted by a three-dimensional (3D) 
granular pack on the container walls and showed that it 
follows an exponential distribution. 2122 Note that, however, 
the coefficient A F does vary with r as shown in Fig. 8. 

3.3 Force chain segment length 

Recent studies have analysed force chain segment lengths 
under pure shear and isotropic compression and found that 
they are exponentially distributed: 23 25 

N cum (l) = Ai exp ’ ( 2 ) 

where N cum (l) is the number of segments of length equal to 
or greater than l. A / and / 0 are fit parameters and / () has di¬ 
mensions of length. The unit of length used is the mean disk 
diameter D = 00l5 2 +001 = 0.0125 m. 

In agreement with previous works, the cumulative number 


distribution for the force chain segment length in this study 
is also found to be exponentially distributed (Fig. 9) with the 
functional form of Eq. (2). The fit parameters are A/ — 9.45 x 
10 2 and Iq = 0.82/1 - 0.01 m at the initial configuration (t - 
0). The characteristic length Iq corresponding to the diameter 
of the small disk is derived from a mere reflection of the effect 
on analysis method. This fact indicates that a segment length 
is meaningless for less than the small disk size. This is natural 
because we consider the force chain structure consisting of 
disks. It is also clearly reflected in Fig. 9 where the steady 
exponential slope is observed only for l > ID. The slope and 
coefficient of the exponential distributions are almost constant 
across all r values, rendering this distribution invariant to the 
manual tapping protocol. 

4. Discussion 

The experimental results discussed thus far show that the 
packing fraction varies with tapping number r through the re¬ 
lation 0 (t) = 0o + A exp(-^), and saturating at an asymptotic 
value of (/>() = 0.79. Additionally, the cumulative distribution 
of force per disk at each r was found to be exponentially dis¬ 
tributed as: N cum {Fd) = A/rexp(-|^). In particular, whereas 
the characteristic force Fq remains invariant to r, the coeffi¬ 
cient A f varies with r (Fig. 8). On the other hand, the cumu¬ 
lative distribution of force chain segment length, which too is 
exponentially distributed (N cum (l) = A/ex p(-j-)), exhibits no 
dependence on tapping number r (Fig. 9). This suggests that 
the evolution of packing fraction 4>(t) leads to the increment 
of internal force within the compacted granular pack, but does 
not lead to creation of new force chain segments. Following 
from these trends, we now explore a speculative relation be¬ 
tween 4>{t ) and the total force F,„, of the granular pack. Since 
(pis a globally averaged structural quantity, it ought to be com- 
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Force chain segment length I [D] 


Fig. 9. Cumulative number distributions of force chain segment lengths in 
log-linear scale. The force chain segment lengths are quoted in units of the 
mean disk diameter, D. The black solid line represents the initial configura¬ 
tion (t = 0) whereas colored lines represent the respective compacted states 
for various r values. The data represent the mean value of nine experimental 
runs. 



Tapping number T 

Fig. 10. Tapping number dependence of the normalized total force F* or 
Mean value of nine runs is reported, and the error bars represent the standard 
error of nine runs. The dashed curve is the fit, F* ot (r ) = F* Q + A* exp(—-^), 
where F* Q = 1.2, A* = -0.2 and r t o = 1.67 are values obtained for the fit 
parameters. 


pared with the total force. 

The total force F to t is defined as F u „ = Ft, where 
Fi is the force per disk on the /th disk, and the summation 
is carried over all disks in the granular pack (k represents 
the total number of disks), with the force threshold placed 
at 0.1 N; forces below this threshold are not included in the 
summation. The total force is measured for the initial con- 



Fig. 11. Relation between F" tot and packing fraction 0. F’ ot is defined by 

F* ot = - — ^ —, where F lot is the sum of force per disk in the granular 

Ftoi(j = 0 ) 

pack. Mean value of nine runs is reported, and the error bars represent the 
standard error of nine runs. The black solid line indicates the linear relation, 
(f; 0 ,(t) - f; o )/(0(t) - <p 0 ) = 14 ± 2. 


figuration and after each manual tap. Accordingly, we de¬ 
fine a normalized total force as F* ot (r) = . where 

F tot (j = 0) represents the total force of the initial configu¬ 
ration. In Fig. 10, we show the normalized total force F* ol as 
a function of r. We can confirm an asymptotic behavior of 
FIJJ) = F *o + A t exp(-; rt) ), where F*„ = 1.2, A* = -0.2, 
and r,o = 1.67. This functional form is similar to that for 
0(r). The comparison of 0 (t) and F* or (r) reveals that t>o — To- 
Therefore, the ratio ( F* ot (r) - F* 0 )/((f>( r) - (f>o ) should be ap¬ 
proximated by A* /A = 14. We independently confirm that 
F* ot (j ) vs 4>(t) scales linearly as shown in Fig. 11. The slope 
of this scaling is 14 + 2 in excellent agreement with the es¬ 
timated result. This linear relation suggests that the process 
of compaction by tapping leads to the increment of granu¬ 
lar internal force, in a linear fashion. This linear dependence 
may possibly result from the optimal amplitude perturbation 
f I -3-4) representing the linear response regime of the sys¬ 
tem. Stronger perturbations may not exhibit a similar depen¬ 
dence between F* oI (t) and This linear relation could be 

potentially useful for estimation of the increment of internal 
force within the compacted granular pack from packing frac¬ 
tion measurements for applications that involve compaction 
processes. 

In this study, we developed the systematic methodology to 
analyse 2D granular pack comprising photoelastic disks. Us¬ 
ing the developed method, granular compaction by manual 
tapping was analysed. Although the interesting structural evo¬ 
lution was revealed in this study, this is still the first step to 
approach granular compaction by tapping using photoelastic 
disks. The result should be compared with the case of con¬ 
trolled tapping using an electromagnetic shaker. Further stud- 
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ies concerning this comparison are in progress at present. The 
result will be published elsewhere. 

5. Conclusion 

In this study, the structural evolution of 2D compacted 
granular pack has been experimentally studied using pho¬ 
toelastic disks. First, we developed the method to measure 
the packing fraction, contact forces, and force chain segment 
lengths using image analysis methods. Then, the dependence 
of these quantities on manual tapping were experimentally 
measured. From the experimental result, the exponentially 
asymptotic behavior of the packing fraction was observed. 
The distributions of applied force per disk and force chain 
segment length at each r were characterized by exponential 
forms. Although the former depends on the tapping number 
r, the latter does not depend on it. The r-dependent total force 
is also shown to exhibit the asymptotic exponential behavior. 
The linear relation of these two functions ((/> and F* ot ) is con¬ 
firmed from the measurements of F* ot and <p. 
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